home *** CD-ROM | disk | FTP | other *** search
/ CD Concept 6 / CD Concept 06.iso / mac / UTILITAIRE / RLaB / rlib / plplot.r < prev    next >
Text File  |  1995-01-07  |  57KB  |  2,535 lines
  1. #
  2. # New plot.r for use with PLPLOT library.
  3. # The help files for these functions are in
  4. # misc/plhelp
  5. #
  6.  
  7. # plplot.r
  8.  
  9. # This file is a part of RLaB ("Our"-LaB)
  10. # Copyright (C) 1994  Ian R. Searle
  11.  
  12. # This program is free software; you can redistribute it and/or modify
  13. # it under the terms of the GNU General Public License as published by
  14. # the Free Software Foundation; either version 2 of the License, or
  15. # (at your option) any later version.
  16.  
  17. # This program is distributed in the hope that it will be useful,
  18. # but WITHOUT ANY WARRANTY; without even the implied warranty of
  19. # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  20. # GNU General Public License for more details.
  21.  
  22. # You should have received a copy of the GNU General Public License
  23. # along with this program; if not, write to the Free Software
  24. # Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
  25.  
  26. # See the file ../COPYING
  27.  
  28. #
  29. # If your system does not deal with Infs and NaNs well, then
  30. # uncomment the following lines.
  31. #
  32.  
  33. # static (isinf, isnan)
  34. #
  35. # isinf = function ( A ) { return (0); };
  36. # isnan = function ( A ) { return (0); };
  37.  
  38. static (WIN)        # The static plot window structure
  39. static (P)        # The active/current plot window
  40.  
  41. #
  42. # Maintain the transformations for 3-D plots.
  43. #
  44.  
  45. static (basex, basey, height)
  46. basex = 2;
  47. basey = 2;
  48. height = 4;
  49.  
  50. #
  51. # Static (private) functions. For use from within
  52. # this file only.
  53. #
  54.  
  55. static (create_plot_object)
  56. static (check_plot_object)
  57. static (xy_scales)
  58. static (x_scales)
  59. static (y_scales)
  60. static (z_scales)
  61. static (XYZ_scales)
  62. static (list_scales)
  63. static (list_sort)
  64. static (hist_scales)
  65. static (plot_matrix)
  66. static (plot_list)
  67. static (check_3d_list)
  68. static (find_char)
  69. static (get_style)
  70. static (make_legend)
  71. static (plhold_first)
  72.  
  73. #
  74. # Defaults
  75. #
  76.  
  77. static (grid_x_default, grid_y_default)
  78. static (grid_3x_default, grid_3y_default, grid_3z_default)
  79.  
  80. grid_x_default = "bcnst";
  81. grid_y_default = "bcnstv";
  82. grid_3x_default = "bnstu";
  83. grid_3y_default = "bnstu";
  84. grid_3z_default = "bcdmnstuv";
  85.  
  86. static (subplot_f)
  87. subplot_f = 0;
  88.  
  89. #
  90. # Create the default plot-object.
  91. # Initialize to all the default values
  92. #
  93.  
  94. if (!exist (WIN)) 
  95. {
  96.   # Create the plot-object list
  97.   WIN = <<>>;
  98. }
  99.  
  100. create_plot_object = function ( N, nx, ny )
  101. {
  102.   if (!exist (N)) { N = 0; }
  103.   
  104.   pobj = <<>>;
  105.   pobj.subplot = 0;        # The current subplot no.
  106.   pobj.nplot = nx*ny;        # Total no. of plots on window
  107.   
  108.   pobj.fontld = 0;        # Loaded extended fonts?
  109.   
  110.   for (i in 1:(nx*ny))
  111.   {
  112.     pobj.style.[i] = "line";    # The type/style of plot to draw
  113.     pobj.pstyle[i] = 1;            # The point-style
  114.     pobj.nbin[i] = 1j;          # The number of bins for a histogram
  115.     pobj.width[i] = 1;        # The pen width for current plot
  116.     pobj.font[i] = 1;        # The current font
  117.     pobj.xlabel[i] = "";
  118.     pobj.ylabel[i] = "";
  119.     pobj.zlabel[i] = "";
  120.     pobj.title[i] = "";
  121.     pobj.orientation[i] = "portrait";
  122.     pobj.desc.[i] = "default";        # The legend description
  123.     pobj.gridx[i] =  grid_x_default;    # Plot axes style, 2D-X
  124.     pobj.gridy[i] =  grid_y_default;    # Plot axes style, 2D-Y
  125.     pobj.grid3x[i] = grid_3x_default;    # Plot axes style, 3D-X
  126.     pobj.grid3y[i] = grid_3y_default;    # Plot axes style, 3D-Y
  127.     pobj.grid3z[i] = grid_3z_default;    # Plot axes style, 3D-Z
  128.     pobj.aspect[i] = 0;                # Plot aspect style
  129.     pobj.alt[i] = 60;
  130.     pobj.az[i] = 45;
  131.     
  132.     pobj.xmin[i] = 1j;
  133.     pobj.xmax[i] = 1j;
  134.     pobj.ymin[i] = 1j;
  135.     pobj.ymax[i] = 1j;
  136.     pobj.zmin[i] = 1j;
  137.     pobj.zmax[i] = 1j;
  138.     
  139.     pobj.page.xp = 0;
  140.     pobj.page.yp = 0;
  141.     pobj.xleng = 400;
  142.     pobj.yleng = 300;
  143.     pobj.xoff = 200;
  144.     pobj.yoff = 200;
  145.  
  146.     for (j in 1:14) { pobj.color[i;j] = j; }       
  147.     for (j in 1:8)  { pobj.lstyle[i;j] = j; }
  148.   }
  149.   
  150.   #
  151.   # Save the newly generated plot-object
  152.   # in a list of plot-objects.
  153.   
  154.   WIN.[N] = pobj;
  155. };
  156.  
  157. ##############################################################################
  158. #
  159. # Check to make sure a plot-object exists. If one
  160. # does not exist, create it.
  161. #
  162.  
  163. check_plot_object = function ()
  164. {
  165.   if (length (WIN) == 0)
  166.   {
  167.     pstart();
  168.     return 0;
  169.   }
  170.   return 1;
  171. };
  172.  
  173. ##############################################################################
  174. #
  175. # Set the current plot window
  176. # Default value = 0
  177. #
  178.  
  179. pwin = function ( N )
  180. {
  181.   check_plot_object ();
  182.   if (!exist (N)) { N = 0; }
  183.   
  184.   # Check to make sure N is valid
  185.   
  186.   for (i in members (WIN))
  187.   {
  188.     if (N == strtod (i))
  189.     {
  190.       _plsstrm (N);
  191.       return P = N;
  192.     }
  193.   }
  194.   printf ("pwin: invalid argument, N = %i\n", N);
  195.   printf ("      valid values are:\n");
  196.   WIN?
  197. };
  198.  
  199. ##############################################################################
  200. #
  201. # Show the current plot-window, and the possibilities
  202. #
  203.  
  204. showpwin = function ( all )
  205. {
  206.   if (length (WIN) == 0)
  207.   {
  208.     printf ("No plot objects\n");
  209.     return 0;
  210.   }
  211.   
  212.   printf ("Current plot-window is:\t\t%i\n", P);
  213.   printf ("Available plot windows are:\t");
  214.   for (i in members (WIN))
  215.   {
  216.     printf ("%s   ", i);
  217.   }
  218.   printf ("\n");
  219.  
  220.   if (exist (all))
  221.   {
  222.     for (i in members (WIN.[P]))
  223.     {
  224.       WIN.[P].[i]
  225.     }
  226.   }
  227. };
  228.  
  229. getplot = function ( win_no )
  230. {
  231.   local (win_no)
  232.  
  233.   if (length (WIN) != 0)
  234.   {
  235.     if (!exist (win_no)) { win_no = P; }
  236.  
  237.     if (exist (WIN.[win_no]))
  238.     {
  239.       return (WIN.[win_no]);
  240.     else
  241.       return 0;
  242.     }
  243.   }
  244.   return <<>>;
  245. };
  246.  
  247.  
  248. ##############################################################################
  249. #
  250. # Set/start/select the plot device
  251. #
  252.  
  253. pstart = function ( nx, ny, dev )
  254. {
  255.   if (!exist (nx)) { nx = 1; }
  256.   if (!exist (ny)) { ny = 1; }
  257.   if (!exist (dev)) { dev = "?"; }
  258.   
  259.   # Create the plot-object
  260.   # First, figure out the index
  261.   if (!exist (P))
  262.   {
  263.     P = 0;
  264.   else
  265.     P = P + 1;
  266.   }
  267.   
  268.   create_plot_object (P, nx, ny);
  269.   _plsstrm (P);
  270.   
  271.   #_plscolbg (255,255,255); # white
  272.   # Default window size for X
  273.   _plspage (0, 0, 400, 300, 200, 200);
  274.   
  275.   # Start up the plot-window
  276.   _plstart (dev, nx, ny);
  277.   
  278.   _plwid (8);
  279.   
  280.   # Turn between plot pause off
  281.   _plspause (0);
  282.   _pltext ();
  283.   
  284.   return P;
  285. };
  286.  
  287. ##############################################################################
  288. #
  289. # Close a plot device. We must destroy the current plot-object
  290. # And switch the output stream back to the default.
  291. #
  292.  
  293. pclose = function ()
  294. {
  295.   if (size (WIN) > 1)
  296.   {   
  297.     #
  298.     # Clear WIN.[P] and reset P to 1st plot-window
  299.     #
  300.  
  301.     clear (WIN.[P]);
  302.     _plend1 ();
  303.     _plsstrm (strtod (members (WIN)[1]));
  304.     P = strtod (members (WIN)[1]);
  305.     return P;
  306.  
  307.   else if (size (WIN) == 1) {
  308.  
  309.     if (exist (WIN.[P])) 
  310.     { 
  311.       clear (WIN.[P]); 
  312.       clear (P);
  313.     }
  314.     _plend1 ();
  315.     return 1;
  316.  
  317.   else if (size (WIN) == 0) {
  318.  
  319.     return 0;
  320.  
  321.   } } }
  322. };
  323.  
  324. ##############################################################################
  325. #
  326. # Close ALL the plot-windows
  327. #
  328.  
  329. pend = function ()
  330. {
  331.   _plend ();
  332.   if (exist (WIN)) { clear (WIN); }
  333.   if (exist (P)) { clear (P); }
  334.   WIN = <<>>;
  335. };
  336.  
  337. ##############################################################################
  338. #
  339. # Change plot aspect ratio
  340. #
  341.  
  342. plaspect = function ( aspect )
  343. {
  344.   check_plot_object ();
  345.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  346.   if (!exist (aspect))
  347.   {
  348.     WIN.[P].aspect[i] = 0;
  349.   else
  350.     WIN.[P].aspect[i] = aspect;
  351.   }
  352. };
  353.  
  354. ##############################################################################
  355. #
  356. # Change plot line style
  357. #
  358.  
  359. plstyle = function ( style )
  360. {
  361.   check_plot_object ();
  362.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  363.  
  364.   if (exist (style)) 
  365.   {
  366.     if (class (style) == "string") 
  367.     {
  368.       WIN.[P].style.[i] = style;
  369.     }
  370.     return 1;
  371.   }
  372.   WIN.[P].style.[i] = "line";  
  373. };
  374.  
  375. ##############################################################################
  376. #
  377. # Get the right value of line-style
  378. #
  379.  
  380. get_style = function ( STY, K )
  381. {
  382.   local (sty);
  383.   sty = mod(K, STY.n);
  384.   if(sty == 0) 
  385.   { 
  386.     sty = STY.n; 
  387.   }
  388.   return STY[sty];
  389. };
  390.  
  391. ##############################################################################
  392. #
  393. # Change fonts
  394. #
  395.  
  396. plfont = function ( font )
  397. {
  398.   check_plot_object ();
  399.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  400.   
  401.   if (!exist (font)) { font = 1; }
  402.   
  403.   if (WIN.[P].fontld == 0)
  404.   {
  405.     _plfontld (1);
  406.     WIN.[P].fontld = 1;
  407.   }
  408.   
  409.   WIN.[P].font[i] = font;
  410.   return P;
  411. };
  412.  
  413. ##############################################################################
  414. #
  415. # Change pen width
  416. #
  417.  
  418. plwid = function ( width )
  419. {
  420.   check_plot_object ();
  421.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  422.   
  423.   if (!exist (width)) { width = 1; }
  424.   WIN.[P].width[i] = width;
  425.   return P;
  426. };
  427.  
  428. ##############################################################################
  429. #
  430. # Place some text on the plot
  431. #
  432.  
  433. plptex = function ( text, x , y , dx , dy , just )
  434. {
  435.   if (!check_plot_object ()) 
  436.   {
  437.     printf ("Must use plot() before plptex()\n");
  438.     return 0;
  439.   }
  440.   
  441.   if (!exist (x)) { x = 0; }
  442.   if (!exist (y)) { y = 0; }
  443.   if (!exist (dx)) { dx = abs(x)+1; }
  444.   if (!exist (dy)) { dy = 0; }
  445.   if (!exist (just)) { just = 0; }
  446.   
  447.   _plptex (x, y, dx, dy, just, text);
  448. };
  449.  
  450. ##############################################################################
  451. #
  452. # Set up the viewing altitude for 3-D plots
  453. #
  454.  
  455. plalt = function ( ALT )
  456. {
  457.   check_plot_object ();
  458.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  459.  
  460.   if (exist (ALT)) 
  461.   { 
  462.     WIN.[P].alt[i] = ALT; 
  463.   else
  464.     WIN.[P].alt[i] = 60;
  465.   }
  466.   return P;
  467. };
  468.  
  469. ##############################################################################
  470. #
  471. # Set the viewing azimuth for 3-D plots
  472. #
  473.  
  474. plaz = function ( AZ )
  475. {
  476.   check_plot_object ();
  477.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  478.  
  479.   if (exist (AZ)) 
  480.   { 
  481.     WIN.[P].az[i] = AZ; 
  482.   else
  483.     WIN.[P].az[i] = 45;
  484.   }
  485.   return P;
  486. };
  487.  
  488. ##############################################################################
  489. #
  490. # Find a character in a string
  491. #
  492.  
  493. find_char = function ( str , char )
  494. {  
  495.   tmp = strsplt (str);
  496.   for (i in 1:tmp.n)
  497.   {
  498.     if (tmp[i] == char) 
  499.     {
  500.       return i;
  501.     }
  502.   }
  503.   return 0;
  504. };
  505.  
  506. ##############################################################################
  507. #
  508. # Sort list element names/labels by numeric order, then string order.
  509. #
  510.  
  511. list_sort = function ( L )
  512. {
  513.   tl = <<>>;
  514.   j = k = 1;
  515.  
  516.   for (i in members (L))
  517.   {
  518.     if (!isnan (strtod (i)))
  519.     {
  520.       num[j] = i;
  521.       j++;
  522.       else
  523.       char[k] = i;
  524.       k++;
  525.     }
  526.   }
  527.   
  528.   # Sort the numeric labels
  529.   
  530.   if (exist (num))
  531.   {
  532.     num = sort (strtod (num)).val;
  533.     tl.num = num;
  534.   }
  535.  
  536.   if (exist (char))
  537.   {
  538.     tl.char = char;
  539.   }
  540.  
  541.   return tl;
  542. };
  543.  
  544. ##############################################################################
  545. #
  546. # Set the subplot, this overides the action in plot().
  547. #
  548.  
  549. subplot = function ( sub )
  550. {
  551.   check_plot_object ();
  552.   
  553.   if (!exist (sub))
  554.   {
  555.     subplot_f = 0;
  556.     _pladv (0);
  557.   else
  558.     if (sub > WIN.[P].nplot)
  559.     {
  560.       error ("Current window does not have this many subplots");
  561.     }
  562.     if (sub > 0)
  563.     {
  564.       WIN.[P].subplot = sub - 1;
  565.       subplot_f = 1;
  566.       _pladv (sub);
  567.     else
  568.       if (sub == 0)
  569.       {
  570.         # Do not advance, stay at current subplot
  571.         WIN.[P].subplot = WIN.[P].subplot - 1;
  572.         subplot_f = 1;
  573.       }
  574.     }
  575.   }
  576. };
  577.  
  578. ##############################################################################
  579. #
  580. # Plot the columns of a matrix (X-Y plot).
  581. #
  582. ##############################################################################
  583.  
  584. plot = function ( data, key )
  585. {
  586.   check_plot_object ();
  587.   
  588.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  589.   if (!exist (key)) { key = 1; }
  590.  
  591.   #
  592.   # Draw the graph
  593.   # Step through the matrix plotting
  594.   # each column versus the 1st
  595.   #
  596.   
  597.   if (class (data) == "num")
  598.   {
  599.     #
  600.     # Set up the plot basics
  601.     #
  602.  
  603.     if (abs (key) > data.nc)
  604.     {
  605.       error_1 ("plot: KEY argument > M.nc");
  606.     }
  607.     
  608.     _plgra ();
  609.     _plcol (1);
  610.     _pllsty (1);
  611.     _plfont (WIN.[P].font[p]);
  612.     _plwid (WIN.[P].width[p]);
  613.  
  614.     if (!subplot_f) 
  615.     {
  616.       _pladv (0);        # Advance 1 subplot
  617.     else
  618.       subplot_f = 0;     # The user has set the subplot
  619.     }
  620.  
  621.     if (WIN.[P].aspect[p] != 0)
  622.     {
  623.       _plvasp (WIN.[P].aspect[p]);
  624.     else
  625.       _plvsta ();
  626.     }
  627.     
  628.     #
  629.     # Compute scale limits
  630.     #
  631.  
  632.     k = find ((1:data.nc) != abs (key));
  633.     if (key > 0)
  634.     {
  635.       if (data.nc != 1)
  636.       {
  637.     x_scales ( real(data)[;key], p, xmin, xmax );
  638.     y_scales ( real(data)[;k],   p, ymin, ymax );
  639.       else
  640.     x_scales ( (1:data.nr)', p, xmin, xmax );
  641.     y_scales ( real(data),   p, ymin, ymax );
  642.       }
  643.     else if (key < 0) {
  644.       x_scales ( real(data)[;k],   p, xmin, xmax );
  645.       y_scales ( real(data)[;abs(key)], p, ymin, ymax );
  646.     else
  647.       x_scales ( (1:data.nr)', p, xmin, xmax );
  648.       y_scales ( real(data),   p, ymin, ymax );
  649.     } }
  650.  
  651.     _plwind (xmin, xmax, ymin, ymax);
  652.     _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  653.     if (plot_matrix ( data, key, p, 0, xmin, xmax, ymin, ymax, ymax-ymin ) < 0) 
  654.     {
  655.       return -1;
  656.     }
  657.     
  658.     else if (class (data) == "list") {
  659.       
  660.       _plgra ();
  661.       _plcol (1);
  662.       _pllsty (1);
  663.       _plfont (WIN.[P].font[p]);
  664.       _plwid (WIN.[P].width[p]);
  665.  
  666.       list_scales ( data, key, p, xmin, xmax, ymin, ymax );
  667.  
  668.       if (!subplot_f) 
  669.       {
  670.     _pladv (0);        # Advance 1 subplot
  671.       else
  672.     subplot_f = 0;     # The user has set the subplot
  673.       }
  674.  
  675.       if (WIN.[P].aspect[p] != 0)
  676.       {
  677.     _plvasp (WIN.[P].aspect[p]);
  678.       else
  679.        _plvsta ();
  680.       }
  681.       
  682.       _plwind (xmin, xmax, ymin, ymax);
  683.       _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  684.       if (plot_list ( data, key, p, xmin, xmax, ymin, ymax ) < 0) 
  685.       { 
  686.     return -1;
  687.       }
  688.  
  689.     else
  690.       error ("plot: un-acceptable argument");
  691.   } }
  692.  
  693.   _pllsty (1);
  694.   _plcol (1);
  695.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  696.   _plflush ();
  697.   _pltext ();
  698.   
  699.   #
  700.   # Increment the plot no. so that next time
  701.   # we use the correct settings.
  702.   #
  703.   
  704.   WIN.[P].subplot = WIN.[P].subplot + 1;
  705.   return P;
  706. };
  707.  
  708. #
  709. # plhold:
  710. # Plot some data, and "hold" on for more.
  711. # Plot the data, setting up the plot as usual the first time.
  712. # On subsequent plots do not do any setup, just plot some
  713. # more. plhold_off must be called to finish up.
  714. #
  715.  
  716. plhold_first = 1;    # True (1) if plhold() has NOT been used.
  717.                         # Or if plhold_off() has been used.
  718.             # False (0) if plhold is in use
  719.  
  720. static (hxmin, hxmax, hymin, hymax)
  721.  
  722. plhold = function ( data, key )
  723. {
  724.   check_plot_object ();
  725.   
  726.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  727.   if (!exist (key)) { key = 1; }
  728.   if (abs (key) > data.nc)
  729.   {
  730.     error_1 ("plot: KEY argument > M.nc");
  731.   }
  732.   
  733.   if (plhold_first)
  734.   {
  735.     if (class (data) == "num")
  736.     {
  737.       #
  738.       # Do the setup ONCE
  739.       #
  740.       hxmin = hxmax = hymin = hymax = 0;
  741.       _plgra ();
  742.       _plcol (1);
  743.       _pllsty (1);
  744.       _plfont (WIN.[P].font[p]);
  745.       _plwid (WIN.[P].width[p]);
  746.  
  747.       if (!subplot_f) 
  748.       {
  749.     _pladv (0);        # Advance 1 subplot
  750.       else
  751.     subplot_f = 0;     # The user has set the subplot
  752.       }
  753.  
  754.       if (WIN.[P].aspect[p] != 0)
  755.       {
  756.     _plvasp (WIN.[P].aspect[p]);
  757.       else
  758.     _plvsta ();
  759.       }
  760.  
  761.       xy_scales ( real(data), p, hxmin, hxmax, hymin, hymax );
  762.       
  763.       k = find ((1:data.nc) != abs (key));
  764.       if (key > 0)
  765.       {
  766.     x_scales ( real(data)[;key], p, xmin, xmax );
  767.     if (data.nc != 1)
  768.     {
  769.       y_scales ( real(data)[;k],   p, ymin, ymax );
  770.         else
  771.       y_scales ( (1:data.nr)',   p, ymin, ymax );
  772.     }
  773.       else if (key < 0) {
  774.     x_scales ( real(data)[;k],   p, xmin, xmax );
  775.     y_scales ( real(data)[;abs(key)], p, ymin, ymax );
  776.       else
  777.     x_scales ( (1:data.nr)', p, xmin, xmax );
  778.     y_scales ( real(data),   p, ymin, ymax );
  779.       } }
  780.  
  781.       _plwind (hxmin, hxmax, hymin, hymax);
  782.       _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  783.       _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);
  784.     else
  785.       error ("plot: un-acceptable argument");
  786.     }
  787.     plhold_first = 0;
  788.   }
  789.   
  790.   if (plot_matrix ( data, key, p, 0, hxmin, hxmax, hymin, hymax, hymax-hymin ) < 0) 
  791.   { 
  792.     return -1; 
  793.   }
  794.  
  795.   _plcol (1);
  796.   _plflush ();
  797.   _pltext ();
  798.  
  799.   return P;
  800. };
  801.  
  802. ##############################################################################
  803. #
  804. # Clean up the plotting environment and get ready
  805. # for normal interactive usage.
  806. #
  807.  
  808. plhold_off = function ( )
  809. {
  810.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  811.   plhold_first = 1;  
  812.   _plcol (1);
  813.   _plflush ();
  814.   _pltext ();
  815.   WIN.[P].subplot = WIN.[P].subplot + 1;
  816.   return P;
  817. };
  818.  
  819. ##############################################################################
  820. #
  821. # Plot a 3-Dimensional graph. The data is composed in a list, with
  822. # elements `x', `y', and `z'. x and y are single-dimension arrays
  823. # (row or column matrices), and z is a two-dimensional array. The
  824. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  825. # the array x can be thought of a "row-labels", and the values of y
  826. # can be thought of as "column-lables" for the 2-dimensioal array z.
  827. #
  828. # At present plot3 can plot 3 distinct lists. Each list may have
  829. # different X, Y, and Z scales.
  830. #
  831. ##############################################################################
  832.  
  833. plot3 = function ( L31, L32, L33 )
  834. {
  835.   check_plot_object ();
  836.  
  837.   #
  838.   # 1st check list contents
  839.   #
  840.   
  841.   if (exist (L31)) { check_3d_list (L31); }
  842.   if (exist (L32)) { check_3d_list (L32); }
  843.   if (exist (L33)) { check_3d_list (L33); }
  844.   
  845.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  846.  
  847.   #
  848.   # Figure out the scale limits. 
  849.   # Needs improvement!
  850.   #
  851.   
  852.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  853.   if (exist (L31)) 
  854.   {
  855.     XYZ_scales (L31.x, L31.y, L31.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  856.     xmin = Xmin; xmax = Xmax;
  857.     ymin = Ymin; ymax = Ymax;
  858.     zmin = Zmin; zmax = Zmax;
  859.   }
  860.   if (exist (L32)) 
  861.   {
  862.     XYZ_scales (L32.x, L32.y, L32.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  863.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  864.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  865.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  866.   }   
  867.   if (exist (L33)) 
  868.   {
  869.     XYZ_scales (L33.x, L33.y, L33.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  870.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  871.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  872.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  873.   }
  874.   
  875.   _plgra ();
  876.   _plcol (1);
  877.   _pllsty (1);
  878.   _plfont (WIN.[P].font[p]);
  879.   _plwid (WIN.[P].width[p]);
  880.   
  881.   # basex = 2; basey = 2; height = 4;
  882.   xmin2d = -2.0; xmax2d = 2.0;
  883.   ymin2d = -3.0; ymax2d = 5.0;
  884.   
  885.   _plenv (xmin2d, xmax2d, ymin2d, ymax2d, 0, -2);
  886.   _plw3d (basex, basey, height, xmin, xmax, ymin, ymax, ...
  887.           zmin, zmax, WIN.[P].alt[p], WIN.[P].az[p]);
  888.   _plbox3 (WIN.[P].grid3x[p], WIN.[P].xlabel[p], 0, 0, ...
  889.            WIN.[P].grid3y[p], WIN.[P].ylabel[p], 0, 0, ...
  890.            WIN.[P].grid3z[p], WIN.[P].zlabel[p], 0, 0);
  891.   _plmtex ("t", 1.0, 0.5, 0.5, WIN.[P].title[p]);
  892.  
  893.   if (exist (L31))
  894.   {
  895.     if (find_char (WIN.[P].grid3x[p], "l"))
  896.     { x = log10 (real (L31.x)); else x = real (L31.x); }
  897.     if (find_char (WIN.[P].grid3y[p], "l"))
  898.     { y = log10 (real (L31.y)); else y = real (L31.y); }
  899.     if (find_char (WIN.[P].grid3z[p], "l"))
  900.     { z = log10 (real (L31.z)); else z = real (L31.z); }
  901.     
  902.     _plcol (2);
  903.     _plot3d (x, y, z, L31.x.n, L31.y.n, 3, 0);
  904.   }
  905.   if (exist (L32))
  906.   {
  907.     if (find_char (WIN.[P].grid3x[p], "l"))
  908.     { x = log10 (real (L32.x)); else x = real (L32.x); }
  909.     if (find_char (WIN.[P].grid3y[p], "l"))
  910.     { y = log10 (real (L32.y)); else y = real (L32.y); }
  911.     if (find_char (WIN.[P].grid3z[p], "l"))
  912.     { z = log10 (real (L32.z)); else z = real (L32.z); }
  913.     
  914.     _plcol (3);
  915.     _pllsty (2);
  916.     _plot3d (x, y, z, L32.x.n, L32.y.n, 3, 0);
  917.   }
  918.   if (exist (L33)) 
  919.   {
  920.     if (find_char (WIN.[P].grid3x[p], "l"))
  921.     { x = log10 (real (L33.x)); else x = real (L33.x); }
  922.     if (find_char (WIN.[P].grid3y[p], "l"))
  923.     { y = log10 (real (L33.y)); else y = real (L33.y); }
  924.     if (find_char (WIN.[P].grid3z[p], "l"))
  925.     { z = log10 (real (L33.z)); else z = real (L33.z); }
  926.     
  927.     _plcol (4);
  928.     _pllsty (3);
  929.     _plot3d (x, y, z, L33.x.n, L33.y.n, 3, 0);
  930.   }
  931.   
  932.   _plflush ();
  933.   _pltext ();
  934.   
  935.   #
  936.   # Increment the plot no. so that next time
  937.   # we use the correct settings.
  938.   #
  939.   
  940.   WIN.[P].subplot = WIN.[P].subplot + 1;
  941.   
  942.   return P;
  943. };
  944.  
  945. ##############################################################################
  946. #
  947. # Plot a 3-Dimensional graph. The data is composed in a list, with
  948. # elements `x', `y', and `z'. x and y are single-dimension arrays
  949. # (row or column matrices), and z is a two-dimensional array. The
  950. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  951. # the array x can be thought of a "row-labels", and the values of y
  952. # can be thought of as "column-lables" for the 2-dimensioal array z.
  953. #
  954. # At present plmesh can plot 3 distinct lists. Each list may have
  955. # different X, Y, and Z scales.
  956. #
  957. ##############################################################################
  958.  
  959. plmesh = function ( L31, L32, L33 )
  960. {
  961.   check_plot_object ();
  962.   
  963.   #
  964.   # 1st check list contents
  965.   #
  966.   
  967.   if (exist (L31)) { check_3d_list (L31); }
  968.   if (exist (L32)) { check_3d_list (L32); }
  969.   if (exist (L33)) { check_3d_list (L33); }
  970.   
  971.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  972.  
  973.   #
  974.   # Figure out the scale limits. 
  975.   # Needs improvement!
  976.   #
  977.   
  978.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  979.   if (exist (L31)) 
  980.   {
  981.     XYZ_scales (L31.x, L31.y, L31.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  982.     xmin = Xmin; xmax = Xmax;
  983.     ymin = Ymin; ymax = Ymax;
  984.     zmin = Zmin; zmax = Zmax;
  985.   }
  986.   if (exist (L32)) 
  987.   {
  988.     XYZ_scales (L32.x, L32.y, L32.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  989.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  990.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  991.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  992.   }   
  993.   if (exist (L33)) 
  994.   {
  995.     XYZ_scales (L33.x, L33.y, L33.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  996.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  997.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  998.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  999.   }
  1000.  
  1001.   _plgra ();
  1002.   _plcol (1);
  1003.   _pllsty (1);
  1004.   _plfont (WIN.[P].font[p]);
  1005.   _plwid (WIN.[P].width[p]);
  1006.   
  1007.   # basex = 2; basey = 2; height = 4;
  1008.   xmin2d = -2.0; xmax2d = 2.0;
  1009.   ymin2d = -3.0; ymax2d = 5.0;
  1010.   
  1011.   _plenv (xmin2d, xmax2d, ymin2d, ymax2d, 0, -2);
  1012.   _plw3d (basex, basey, height, xmin, xmax, ymin, ymax, ...
  1013.           zmin, zmax, WIN.[P].alt[p], WIN.[P].az[p]);
  1014.   _plbox3 (WIN.[P].grid3x[p], WIN.[P].xlabel[p], 0, 0, ...
  1015.            WIN.[P].grid3y[p], WIN.[P].ylabel[p], 0, 0, ...
  1016.            WIN.[P].grid3z[p], WIN.[P].zlabel[p], 0, 0);
  1017.   _plmtex ("t", 1.0, 0.5, 0.5, WIN.[P].title[p]);
  1018.  
  1019.   if (exist (L31))
  1020.   {
  1021.     if (find_char (WIN.[P].grid3x[p], "l"))
  1022.     { x = log10 (real (L31.x)); else x = real (L31.x); }
  1023.     if (find_char (WIN.[P].grid3y[p], "l"))
  1024.     { y = log10 (real (L31.y)); else y = real (L31.y); }
  1025.     if (find_char (WIN.[P].grid3z[p], "l"))
  1026.     { z = log10 (real (L31.z)); else z = real (L31.z); }
  1027.     
  1028.     _plcol (2);
  1029.     _plmesh (x, y, z, L31.x.n, L31.y.n, 3);
  1030.   }
  1031.   if (exist (L32))
  1032.   {
  1033.     if (find_char (WIN.[P].grid3x[p], "l"))
  1034.     { x = log10 (real (L32.x)); else x = real (L32.x); }
  1035.     if (find_char (WIN.[P].grid3y[p], "l"))
  1036.     { y = log10 (real (L32.y)); else y = real (L32.y); }
  1037.     if (find_char (WIN.[P].grid3z[p], "l"))
  1038.     { z = log10 (real (L32.z)); else z = real (L32.z); }
  1039.     
  1040.     _plcol (3);
  1041.     _pllsty (2);
  1042.     _plmesh (x, y, z, L32.x.n, L32.y.n, 3);
  1043.   }
  1044.   if (exist (L33)) 
  1045.   {
  1046.     if (find_char (WIN.[P].grid3x[p], "l"))
  1047.     { x = log10 (real (L33.x)); else x = real (L33.x); }
  1048.     if (find_char (WIN.[P].grid3y[p], "l"))
  1049.     { y = log10 (real (L33.y)); else y = real (L33.y); }
  1050.     if (find_char (WIN.[P].grid3z[p], "l"))
  1051.     { z = log10 (real (L33.z)); else z = real (L33.z); }
  1052.     
  1053.     _plcol (4);
  1054.     _pllsty (3);
  1055.     _plmesh (x, y, z, L33.x.n, L33.y.n, 3);
  1056.   }
  1057.  
  1058.   _plflush ();
  1059.   _pltext ();
  1060.   
  1061.   #
  1062.   # Increment the plot no. so that next time
  1063.   # we use the correct settings.
  1064.   #
  1065.   
  1066.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1067.   
  1068.   return P;
  1069. };
  1070.  
  1071. ##############################################################################
  1072. #
  1073. # Plot contours. The data is composed in a list, with
  1074. # elements `x', `y', and `z'. x and y are single-dimension arrays
  1075. # (row or column matrices), and z is a two-dimensional array. The
  1076. # array z, is a function of x and y: z = f(x,y). Thus, the values in
  1077. # the array x can be thought of a "row-labels", and the values of y
  1078. # can be thought of as "column-lables" for the 2-dimensioal array z.
  1079. #
  1080. ##############################################################################
  1081.  
  1082. plcont = function ( CL )
  1083. {
  1084.   check_plot_object ();
  1085.   
  1086.   #
  1087.   # 1st check list contents
  1088.   #
  1089.   
  1090.   if (exist (CL)) { check_3d_list (CL); }
  1091.  
  1092.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1093.   
  1094.   #
  1095.   # Figure out the scale limits. 
  1096.   # Needs improvement!
  1097.   #
  1098.  
  1099.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  1100.   if (exist (CL)) 
  1101.   {
  1102.     XYZ_scales (CL.x, CL.y, CL.z, p, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax);
  1103.     if (Xmin < xmin) { xmin = Xmin; } if (Xmax > xmax) { xmax = Xmax; }
  1104.     if (Ymin < ymin) { ymin = Ymin; } if (Ymax > ymax) { ymax = Ymax; }
  1105.     if (Zmin < zmin) { zmin = Zmin; } if (Zmax > zmax) { zmax = Zmax; }
  1106.   }
  1107.  
  1108.   _plgra ();
  1109.   _plcol (1);
  1110.   _pllsty (1);
  1111.   _plfont (WIN.[P].font[p]);
  1112.   _plwid (WIN.[P].width[p]);
  1113.  
  1114.   #
  1115.   # Set up the 1st viewport for drawing the plot.
  1116.   #
  1117.  
  1118.   if (!subplot_f) 
  1119.   {
  1120.     _pladv (0);        # Advance 1 subplot
  1121.   else
  1122.     subplot_f = 0;     # The user has set the subplot
  1123.   }
  1124.  
  1125.   _plvpas (0.15, 0.75, 0.15, 0.85, WIN.[P].aspect[p]);
  1126.   # _plvpor (0.15, 0.75, 0.15, 0.85);
  1127.   _plwind (xmin, xmax, ymin, ymax);
  1128.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  1129.  
  1130.   # Convert the data to log data if necessary.
  1131.   if (find_char (WIN.[P].gridx[p], "l"))
  1132.   { x = log10 (real (CL.x)); else x = real (CL.x); }
  1133.   if (find_char (WIN.[P].gridy[p], "l"))
  1134.   { y = log10 (real (CL.y)); else y = real (CL.y); }
  1135.   z = real (CL.z);
  1136.   
  1137.   if (exist (CL.clevel))
  1138.   {
  1139.     clevel = CL.clevel;
  1140.   else
  1141.     clevel = linspace(zmin, zmax, 10);
  1142.   }
  1143.  
  1144.   #
  1145.   # Draw the contours
  1146.   #
  1147.  
  1148.   l = 1;
  1149.   for (i in 1:clevel.n)
  1150.   {
  1151.     k = mod (i-1, 14) + 1;
  1152.     j = mod (i-1, 8) + 1;
  1153.     _pllsty(j);
  1154.     _plcol (1+k);
  1155.     if (_plcont (x, y, z, 1, CL.x.n, 1, CL.y.n, clevel[i]))
  1156.     {
  1157.       llevel[l] = clevel[i];
  1158.       l = l + 1;
  1159.     }
  1160.   }
  1161.  
  1162.   #
  1163.   # Reset color and draw the labels.
  1164.   #
  1165.  
  1166.   _plcol (1);
  1167.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  1168.  
  1169.   #
  1170.   # Draw the contour legend. Use a new viewport to the right
  1171.   # of the contour plot.
  1172.   #
  1173.  
  1174.   #_plvpas (0.75, 1.0, 0.15, 0.85, WIN.[P].aspect[p]);
  1175.   _plvpor (0.75, 1.0, 0.15, 0.85);
  1176.   _plwind (0, 1, 0, 1);
  1177.  
  1178.   v = 1 - 1/(2*llevel.n);
  1179.  
  1180.   for (i in 1:llevel.n)
  1181.   {
  1182.     xl = [0.1, 0.2, 0.3]';
  1183.     yl = [v, v, v]';
  1184.     v = v - 1/llevel.n;
  1185.  
  1186.     k = mod (i-1, 14) + 1;
  1187.     j = mod (i-1, 8) + 1;
  1188.  
  1189.     _plcol (1+k);
  1190.     _pllsty (j);
  1191.  
  1192.     _plline (3, xl, yl);
  1193.     sprintf (stmp, "%.2g", llevel[i]);
  1194.     plptex (stmp, xl[3]+.1, yl[3], , , 0);
  1195.   }
  1196.  
  1197.   # Flush  and go back to text mode.
  1198.   _plflush ();
  1199.   _pltext ();
  1200.   
  1201.   #
  1202.   # Increment the plot no. so that next time
  1203.   # we use the correct settings.
  1204.   #
  1205.   
  1206.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1207.   
  1208.   return P;
  1209. };
  1210.  
  1211. ##############################################################################
  1212. #
  1213. # Plot 3-D lines, etc...
  1214. #
  1215.  
  1216. pl3d = function ( X, Y, Z, BR )
  1217. {
  1218.   local (X, Y, Z, BR)
  1219.   check_plot_object ();
  1220.   
  1221.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1222.  
  1223.   #
  1224.   # Some basic checks
  1225.   #
  1226.  
  1227.   if ((N = X.n) != Y.n) { error ("pl3d: X and Y must have same length"); }
  1228.   if (N != Z.n) { error ("pl3d: X and Z must have same length"); }
  1229.  
  1230.   if (!exist (BR)) { BR = N; }
  1231.   if (mod (N, BR) != 0) { error ("pl3d: X.n must be divisible by BR"); }
  1232.   iBR = int (N / BR);
  1233.   if (iBR == 1) { k = N; else k = BR; }
  1234.  
  1235.   #
  1236.   # Figure out the scale limits. 
  1237.   # Needs improvement!
  1238.   #
  1239.   
  1240.   xmin = xmax = ymin = ymax = zmin = zmax = 0;
  1241.   XYZ_scales (X, Y, Z, p, xmin, xmax, ymin, ymax, zmin, zmax);
  1242.  
  1243.   _plgra ();
  1244.   _plcol (1);
  1245.   _pllsty (1);
  1246.   _plfont (WIN.[P].font[p]);
  1247.   _plwid (WIN.[P].width[p]);
  1248.   
  1249.   if (!subplot_f) 
  1250.   {
  1251.     _pladv (0);        # Advance 1 subplot
  1252.     else
  1253.     subplot_f = 0;     # The user has set the subplot
  1254.   }
  1255.  
  1256.   if (find_char (WIN.[P].grid3x[p], "l"))
  1257.   { X = log10 (real (X)); else X = real (X); }
  1258.   if (find_char (WIN.[P].grid3y[p], "l"))
  1259.   { Y = log10 (real (Y)); else Y = real (Y); }
  1260.   if (find_char (WIN.[P].grid3z[p], "l"))
  1261.   { Z = log10 (real (Z)); else Z = real (Z); }
  1262.     
  1263.   # basex = 2; basey = 2; height = 4;
  1264.   xmin2d = -2.0; xmax2d = 2.0;
  1265.   ymin2d = -3.0; ymax2d = 5.0;
  1266.   
  1267.   _plvasp (WIN.[P].aspect[p]);
  1268.   _plwind (xmin2d, xmax2d, ymin2d, ymax2d);
  1269.   _plw3d (basex, basey, height, xmin, xmax, ymin, ymax, ...
  1270.           zmin, zmax, WIN.[P].alt[p], WIN.[P].az[p]);
  1271.   _plbox3 (WIN.[P].grid3x[p], WIN.[P].xlabel[p], 0, 0, ...
  1272.            WIN.[P].grid3y[p], WIN.[P].ylabel[p], 0, 0, ...
  1273.            WIN.[P].grid3z[p], WIN.[P].zlabel[p], 0, 0);
  1274.   _plmtex ("t", 1.0, 0.5, 0.5, WIN.[P].title[p]);
  1275.  
  1276.   _plcol (2);
  1277.   for (i in 1:iBR)
  1278.   {
  1279.     j = [(i-1)*k+1:i*k];
  1280.     _plline3 (k, X[j], Y[j], Z[j]);
  1281.   }
  1282.   _plflush ();
  1283.   _pltext ();
  1284.   
  1285.   #
  1286.   # Increment the plot no. so that next time
  1287.   # we use the correct settings.
  1288.   #
  1289.   
  1290.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1291.   
  1292.   return P;
  1293. };
  1294.  
  1295. ##############################################################################
  1296. #
  1297. # error bar plot
  1298. #
  1299.  
  1300. plerry = function (x, y, y_low, y_high)
  1301. {
  1302.   check_plot_object ();
  1303.   
  1304.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1305.   WIN.[P].desc.[p] = 1j;    // use own legend
  1306.   
  1307.   if (x.nr != y.nr || x.nr != y_low.nr || x.nr != y_high.nr) 
  1308.   {
  1309.     error(" Size inconsistent in plerry.");
  1310.   }
  1311.   _plgra ();
  1312.   _plcol (1);
  1313.   _pllsty (1);
  1314.   _plfont (WIN.[P].font[p]);
  1315.   _plwid (WIN.[P].width[p]);
  1316.   xy_scales ( real([x,y,y_low,y_high]), p, xmin, xmax, ymin, ymax );
  1317.   
  1318.   if (!subplot_f) 
  1319.   {
  1320.     _pladv (0);        # Advance 1 subplot
  1321.   else
  1322.     subplot_f = 0;     # The user has set the subplot
  1323.   }
  1324.  
  1325.   if (WIN.[P].aspect[p] != 0)
  1326.   {
  1327.     _plvasp (WIN.[P].aspect[p]);
  1328.   else
  1329.     _plvsta ();
  1330.   }
  1331.  
  1332.   _plwind (xmin, xmax, ymin, ymax);
  1333.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  1334.  
  1335.   if (plot_matrix ( [x,y], p, 0, xmin, xmax, ymin, ymax, ymax-ymin ) < 0) 
  1336.   { 
  1337.     return -1; 
  1338.   }
  1339.  
  1340.   _plcol (3);
  1341.   _plerry(x.nr, x, y_low, y_high);
  1342.   _plcol (1);
  1343.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  1344.   _plflush ();
  1345.   _pltext ();
  1346.   
  1347.   #
  1348.   # Increment the plot no. so that next time
  1349.   # we use the correct settings.
  1350.   #
  1351.   
  1352.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1353.   return P;      
  1354. };
  1355.  
  1356. ##############################################################################
  1357. #
  1358. # Plot a Histogram(s), from the columns of a matrix.
  1359. #
  1360. ##############################################################################
  1361.  
  1362. plhist = function ( M , nbin )
  1363. {
  1364.   check_plot_object ();
  1365.   
  1366.   if (!exist (nbin)) { nbin = 10; }
  1367.   
  1368.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1369.   np = M.nr;
  1370.   
  1371.   # Compute max/min values of data
  1372.   
  1373.   ymin = min (min (real (M)));
  1374.   ymax = max (max (real (M)));
  1375.   
  1376.   #
  1377.   # Check computed scale limits against user's
  1378.   #
  1379.   
  1380.   if (WIN.[P].ymin[p] != 1j) { ymin = WIN.[P].ymin[p]; }
  1381.   if (WIN.[P].ymax[p] != 1j) { ymax = WIN.[P].ymax[p]; }
  1382.  
  1383.   _plgra ();
  1384.   _plcol (1);
  1385.   _plfont (WIN.[P].font[p]);
  1386.   _plwid (WIN.[P].width[p]);
  1387.   
  1388.   for (i in 1:M.nc) 
  1389.   { 
  1390.     hscale[i] = hist_scales (M[;i], nbin);
  1391.   }  
  1392.  
  1393.   if (!subplot_f) {
  1394.     _pladv (0);        # Advance 1 subplot
  1395.   else
  1396.     subplot_f = 0;     # The user has set the subplot
  1397.   }
  1398.  
  1399.   if (WIN.[P].aspect[p] != 0)
  1400.   {
  1401.     _plvasp (WIN.[P].aspect[p]);
  1402.   else
  1403.     _plvsta ();
  1404.   }
  1405.  
  1406.   _plwind (ymin, ymax, 0, max (hscale));
  1407.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  1408.  
  1409.   v = max (hscale);
  1410.   xmax = ymax;
  1411.   for (i in 1:M.nc)
  1412.   {
  1413.     k = mod (i, 14) + 1;
  1414.     _plcol (WIN.[P].color[p;k]);
  1415.     _plhist (np, real(M[;i]), ymin, ymax, nbin, 1);
  1416.     
  1417.     if (!any (any (WIN.[P].desc.[p] == 1j)))
  1418.     {
  1419.       # Use the default if necessary
  1420.       if (WIN.[P].desc.[p][1] == "default") 
  1421.       {
  1422.     desc = "c"+num2str(i);
  1423.       else if (WIN.[P].desc.[p].n >= i) {
  1424.     desc = WIN.[P].desc.[p][i];
  1425.       else
  1426.     # Not sure what to do, user has messed up.
  1427.     desc = "";
  1428.       } }
  1429.  
  1430.       v = v - max(hscale)/11;
  1431.       xl = (ymax-ymin)*[10.5/12, 11/12, 11.5/12]' + ymin;
  1432.       yl = [v, v, v]';
  1433.  
  1434.       _plline (3, xl, yl);
  1435.       plptex(desc, xl[1]-(ymax-ymin)/25, yl[3], , , 1);
  1436.     }
  1437.   }
  1438.  
  1439.   _plcol (1);
  1440.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  1441.   _plflush ();
  1442.   _pltext ();
  1443.  
  1444.   #
  1445.   # Increment the plot no. so that next time
  1446.   # we use the correct settings.
  1447.   #
  1448.   
  1449.   WIN.[P].subplot = WIN.[P].subplot + 1;
  1450.  
  1451.   return 1;
  1452. };
  1453.  
  1454. ##############################################################################
  1455. #
  1456. # Various support functions for the WIN list
  1457. #
  1458. ##############################################################################
  1459.  
  1460. #
  1461. # Replot
  1462. #
  1463.  
  1464. replot = function ( )
  1465. {
  1466.   check_plot_object ();
  1467.   _replot ();
  1468. };
  1469.  
  1470. ##############################################################################
  1471. #
  1472. # Set the X-axis label
  1473. #
  1474.  
  1475. xlabel = function ( xstr )
  1476. {
  1477.   check_plot_object ();
  1478.   if (!exist (xstr)) { xstr = ""; }
  1479.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1480.   WIN.[P].xlabel[i] = xstr;
  1481. };
  1482.  
  1483. ##############################################################################
  1484. #
  1485. # Set the Y-axis label
  1486. #
  1487.  
  1488. ylabel = function ( xstr )
  1489. {
  1490.   check_plot_object ();
  1491.   if (!exist (xstr)) { xstr = ""; }
  1492.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1493.   WIN.[P].ylabel[i] = xstr;
  1494. };
  1495.  
  1496. ##############################################################################
  1497. #
  1498. # Set the Z-axis label
  1499. #
  1500.  
  1501. zlabel = function ( xstr )
  1502. {
  1503.   check_plot_object ();
  1504.   if (!exist (xstr)) { xstr = ""; }
  1505.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1506.   WIN.[P].zlabel[i] = xstr;
  1507. };
  1508.  
  1509. ##############################################################################
  1510. #
  1511. # Set the plot-title
  1512. #
  1513.  
  1514. ptitle = function ( xstr )
  1515. {
  1516.   check_plot_object ();
  1517.   if (!exist (xstr)) { xstr = ""; }
  1518.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  1519.   WIN.[P].title[p] = xstr;
  1520. };
  1521.  
  1522. ##############################################################################
  1523. #
  1524. # Set the scale limits.
  1525. #
  1526.  
  1527. plimits = function ( xmin, xmax, ymin, ymax, zmin, zmax )
  1528. {
  1529.   check_plot_object ();
  1530.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1531.  
  1532.   if (exist (xmin)) 
  1533.   {
  1534.     WIN.[P].xmin[i] = xmin;
  1535.   else
  1536.     WIN.[P].xmin[i] = 1j;
  1537.   }
  1538.   if (exist (xmax)) 
  1539.   { 
  1540.     WIN.[P].xmax[i] = xmax;
  1541.   else
  1542.     WIN.[P].xmax[i] = 1j;
  1543.   }
  1544.  
  1545.   if (exist (ymin)) 
  1546.   {
  1547.     WIN.[P].ymin[i] = ymin;
  1548.   else
  1549.     WIN.[P].ymin[i] = 1j;
  1550.   }
  1551.   if (exist (ymax)) 
  1552.   {
  1553.     WIN.[P].ymax[i] = ymax;
  1554.   else
  1555.     WIN.[P].ymax[i] = 1j;
  1556.   }
  1557.  
  1558.   if (exist (zmin)) 
  1559.   {
  1560.     WIN.[P].zmin[i] = zmin;
  1561.   else
  1562.     WIN.[P].zmin[i] = 1j;
  1563.   }
  1564.   if (exist (zmax)) 
  1565.   {
  1566.     WIN.[P].zmax[i] = zmax;
  1567.   else
  1568.     WIN.[P].zmax[i] = 1j;
  1569.   }
  1570. };
  1571.  
  1572. ##############################################################################
  1573. #
  1574. # Set 2-D grid styles. A not-so-friendly interface.
  1575. #
  1576.  
  1577. plgrid = function ( sty_x, sty_y )
  1578. {
  1579.   check_plot_object ();
  1580.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1581.  
  1582.   if (exist (sty_x)) 
  1583.   { 
  1584.     if (class (sty_x) == "string")
  1585.     {
  1586.       WIN.[P].gridx[i] = sty_x;
  1587.     else
  1588.       error ("plgrid: requires string argument GRID_STY_X");
  1589.     }
  1590.   else
  1591.     WIN.[P].gridx[i] = grid_x_default;
  1592.   }
  1593.   if (exist (sty_y)) 
  1594.   { 
  1595.     if (class (sty_y) == "string")
  1596.     {
  1597.       WIN.[P].gridy[i] = sty_y;
  1598.     else
  1599.       error ("plgrid: requires string argument GRID_STY_Y");
  1600.     }
  1601.   else
  1602.     WIN.[P].gridy[i] = grid_y_default;
  1603.   }
  1604. };
  1605.  
  1606. ##############################################################################
  1607. #
  1608. # Set 3-D grid (axis) styles
  1609. #
  1610.  
  1611. plgrid3 = function ( sty_x, sty_y, sty_z )
  1612. {  
  1613.   check_plot_object ();
  1614.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1615.   if (exist (sty_x)) 
  1616.   { 
  1617.     if (class (sty_x) == "string")
  1618.     {
  1619.       WIN.[P].grid3x[i] = sty_x;
  1620.     else
  1621.       error ("plgrid3: requires string argument GRID_STY_X");
  1622.     }
  1623.   else
  1624.     WIN.[P].grid3x[i] = grid_3x_default;
  1625.   }
  1626.   if (exist (sty_y)) 
  1627.   { 
  1628.     if (class (sty_y) == "string")
  1629.     {
  1630.       WIN.[P].grid3y[i] = sty_y;
  1631.     else
  1632.       error ("plgrid3: requires string argument GRID_STY_Y");
  1633.     }
  1634.   else
  1635.     WIN.[P].grid3y[i] = grid_3y_default;
  1636.   }
  1637.   if (exist (sty_z)) 
  1638.   { 
  1639.     if (class (sty_z) == "string")
  1640.     {
  1641.       WIN.[P].grid3z[i] = sty_z;
  1642.     else
  1643.       error ("plgrid3: requires string argument GRID_STY_Z");
  1644.     }
  1645.   else
  1646.     WIN.[P].grid3z[i] = grid_3z_default;
  1647.   }
  1648. };
  1649.  
  1650. ##############################################################################
  1651. #
  1652. # A friendlier interface to changing 2-D grid/axis
  1653. # styles.
  1654. #
  1655.  
  1656. plaxis = function ( X_STR, Y_STR )
  1657. {
  1658.   check_plot_object ();
  1659.   i = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;
  1660.   
  1661.   if (exist (X_STR))
  1662.   {
  1663.     if (X_STR == "log") { WIN.[P].gridx[i] = "bcngstl"; }
  1664.   else
  1665.     WIN.[P].gridx[i] = grid_x_default;
  1666.   }
  1667.   
  1668.   if (exist (Y_STR))
  1669.   {
  1670.     if (Y_STR == "log") { WIN.[P].gridy[i] = "bcngstlv"; }
  1671.   else
  1672.     WIN.[P].gridy[i] = grid_y_default;
  1673.   }
  1674.   return P;
  1675. };
  1676.  
  1677. ##############################################################################
  1678. #
  1679. # Various internal support functions. Eventually these will be static.
  1680. #
  1681. ##############################################################################
  1682.  
  1683. #
  1684. # Find the X or Y scale limits .
  1685. # M can be a multi-column matrix, all columns
  1686. # will be used.
  1687. #
  1688.  
  1689. x_scales = function ( M, p, xmin, xmax )
  1690. {  
  1691.   #
  1692.   # 1st check for un-plottable data
  1693.   #
  1694.   
  1695.   if (any (any (isinf (M))))
  1696.   { error ("plot: cannot plot Infs"); }
  1697.   if (any (any (isnan (M))))
  1698.   { error ("plot: cannot plot NaNs"); }
  1699.   
  1700.   xmin = min (min (M));
  1701.   xmax = max (max (M));
  1702.   
  1703.   #
  1704.   # Check computed scale limits against user's
  1705.   #
  1706.   
  1707.   if (WIN.[P].xmin[p] != 1j) { xmin = WIN.[P].xmin[p]; }
  1708.   if (WIN.[P].xmax[p] != 1j) { xmax = WIN.[P].xmax[p]; }
  1709.  
  1710.   #
  1711.   # Check for potential errors
  1712.   #
  1713.   
  1714.   if (xmin == xmax) 
  1715.   { 
  1716.     # As good a guess as any
  1717.     xmin = xmin - 1;
  1718.     xmax = xmax + 1;
  1719.   }
  1720.   
  1721.   #
  1722.   # Finally, adjust if log-scales
  1723.   #
  1724.   
  1725.   if (find_char (WIN.[P].gridx[p], "l"))
  1726.   {
  1727.     if (xmin <= 0 || xmax <= 0) { error ("cannot plot log(x<=0)"); }
  1728.     xmin = log10 (xmin);
  1729.     xmax = log10 (xmax);
  1730.   }
  1731.   
  1732.   return 1;
  1733. };
  1734.  
  1735. y_scales = function ( M, p, xmin, xmax )
  1736. {
  1737.   
  1738.   #
  1739.   # 1st check for un-plottable data
  1740.   #
  1741.   
  1742.   if (any (any (isinf (M))))
  1743.   { error ("plot: cannot plot Infs"); }
  1744.   if (any (any (isnan (M))))
  1745.   { error ("plot: cannot plot NaNs"); }
  1746.   
  1747.   xmin = min (min (M));
  1748.   xmax = max (max (M));
  1749.   
  1750.   #
  1751.   # Check computed scale limits against user's
  1752.   #
  1753.   
  1754.   if (WIN.[P].ymin[p] != 1j) { xmin = WIN.[P].ymin[p]; }
  1755.   if (WIN.[P].ymax[p] != 1j) { xmax = WIN.[P].ymax[p]; }
  1756.  
  1757.   #
  1758.   # Check for potential errors
  1759.   #
  1760.   
  1761.   if (xmin == xmax) 
  1762.   { 
  1763.     # As good a guess as any
  1764.     xmin = xmin - 1;
  1765.     xmax = xmax + 1;
  1766.   }
  1767.   
  1768.   #
  1769.   # Finally, adjust if log-scales
  1770.   #
  1771.   
  1772.   if (find_char (WIN.[P].gridy[p], "l"))
  1773.   {
  1774.     if (xmin <= 0 || xmax <= 0) { error ("cannot plot log(y<=0)"); }
  1775.     xmin = log10 (xmin);
  1776.     xmax = log10 (xmax);
  1777.   }
  1778.   
  1779.   return 1;
  1780. };
  1781.  
  1782. z_scales = function ( M, p, xmin, xmax )
  1783. {
  1784.   
  1785.   #
  1786.   # 1st check for un-plottable data
  1787.   #
  1788.   
  1789.   if (any (any (isinf (M))))
  1790.   { error ("plot: cannot plot Infs"); }
  1791.   if (any (any (isnan (M))))
  1792.   { error ("plot: cannot plot NaNs"); }
  1793.   
  1794.   xmin = min (min (M));
  1795.   xmax = max (max (M));
  1796.   
  1797.   #
  1798.   # Check computed scale limits against user's
  1799.   #
  1800.   
  1801.   if (WIN.[P].zmin[p] != 1j) { xmin = WIN.[P].zmin[p]; }
  1802.   if (WIN.[P].zmax[p] != 1j) { xmax = WIN.[P].zmax[p]; }
  1803.  
  1804.   #
  1805.   # Check for potential errors
  1806.   #
  1807.   
  1808.   if (xmin == xmax) 
  1809.   { 
  1810.     # As good a guess as any
  1811.     xmin = xmin - 1;
  1812.     xmax = xmax + 1;
  1813.   }
  1814.   
  1815.   #
  1816.   # Finally, adjust if log-scales
  1817.   #
  1818.   
  1819.   if (find_char (WIN.[P].gridz[p], "l"))
  1820.   {
  1821.     if (xmin <= 0 || xmax <= 0) { error ("plot: cannot plot log(z<=0)"); }
  1822.     xmin = log10 (xmin);
  1823.     xmax = log10 (xmax);
  1824.   }
  1825.   
  1826.   return 1;
  1827. };
  1828.  
  1829. ##############################################################################
  1830. #
  1831. # Find the X and Y scales for a single matrix. (OLD)
  1832. #
  1833.  
  1834. xy_scales = function ( M, p, xmin, xmax, ymin, ymax )
  1835. {  
  1836.   #
  1837.   # 1st check for un-plottable data
  1838.   #
  1839.   
  1840.   if (any (any (isinf (M))))
  1841.   { error ("plot: cannot plot infs"); }
  1842.   if (any (any (isnan (M))))
  1843.   { error ("plot: cannot plot NaNs"); }
  1844.   
  1845.   if (M.nc == 1)
  1846.   {
  1847.     xmin = 1;
  1848.     xmax = M.nr;
  1849.     ymin = min (M);
  1850.     ymax = max (M);
  1851.   else
  1852.     xmin = min (M[;1]);
  1853.     xmax = max (M[;1]);
  1854.     ymin = min (min (M[;2:M.nc]));
  1855.     ymax = max (max (M[;2:M.nc]));
  1856.   }
  1857.   
  1858.   #
  1859.   # Check computed scale limits against user's
  1860.   #
  1861.   
  1862.   if (WIN.[P].xmin[p] != 1j) { xmin = WIN.[P].xmin[p]; }
  1863.   if (WIN.[P].xmax[p] != 1j) { xmax = WIN.[P].xmax[p]; }
  1864.   if (WIN.[P].ymin[p] != 1j) { ymin = WIN.[P].ymin[p]; }
  1865.   if (WIN.[P].ymax[p] != 1j) { ymax = WIN.[P].ymax[p]; }
  1866.   
  1867.   #
  1868.   # Check for potential errors
  1869.   #
  1870.   
  1871.   if (xmin == xmax) 
  1872.   { 
  1873.     xmin = xmin - 1;
  1874.     xmax = xmax + 1;
  1875.   }
  1876.   
  1877.   if (ymin == ymax)
  1878.   {
  1879.     ymin = ymin - 1;
  1880.     ymax = ymax + 1;
  1881.   }
  1882.   
  1883.   #
  1884.   # Finally, adjust if log-scales
  1885.   #
  1886.   
  1887.   if (find_char (WIN.[P].gridx[p], "l"))
  1888.   {
  1889.     if (xmin <= 0 || xmax <= 0) { error ("plot: cannot plot log <= 0"); }
  1890.     xmin = log10 (xmin);
  1891.     xmax = log10 (xmax);
  1892.   }
  1893.   if (find_char (WIN.[P].gridy[p], "l"))
  1894.   {
  1895.     if (ymin <= 0 || ymax <= 0) { error ("plot: cannot plot log <= 0"); }
  1896.     ymin = log10 (ymin);
  1897.     ymax = log10 (ymax);
  1898.   }
  1899.   
  1900.   return 1;
  1901. };
  1902.  
  1903. ##############################################################################
  1904. #
  1905. # Find the X, Y and Z scales for a single matrix.
  1906. #
  1907.  
  1908. XYZ_scales = function ( X, Y, Z, p, xmin, xmax, ymin, ymax, zmin, zmax )
  1909. {
  1910.   # X - scale
  1911.   if (any (any (isinf (X))))
  1912.   { error ("cannot plot infs"); }
  1913.   if (any (any (isnan (X))))
  1914.   { error ("cannot plot NaNs"); }
  1915.   
  1916.   xmin = min (real (X));
  1917.   xmax = max (real (X));
  1918.   
  1919.   # Y - scale
  1920.   if (any (any (isinf (Y))))
  1921.   { error ("cannot plot infs"); }
  1922.   if (any (any (isnan (Y))))
  1923.   { error ("cannot plot NaNs"); }
  1924.   
  1925.   ymin = min (real (Y));
  1926.   ymax = max (real (Y));
  1927.   
  1928.   # Z - scale
  1929.   if (any (any (isinf (Y))))
  1930.   { error ("cannot plot infs"); }
  1931.   if (any (any (isnan (Y))))
  1932.   { error ("cannot plot NaNs"); }
  1933.   
  1934.   zmin = min (min (real (Z)));
  1935.   zmax = max (max (real (Z)));
  1936.   
  1937.   #
  1938.   # Check computed scale limits against user's
  1939.   #
  1940.   
  1941.   if (WIN.[P].xmin[p] != 1j) { xmin = WIN.[P].xmin[p]; }
  1942.   if (WIN.[P].xmax[p] != 1j) { xmax = WIN.[P].xmax[p]; }
  1943.   if (WIN.[P].ymin[p] != 1j) { ymin = WIN.[P].ymin[p]; }
  1944.   if (WIN.[P].ymax[p] != 1j) { ymax = WIN.[P].ymax[p]; }
  1945.   if (WIN.[P].zmin[p] != 1j) { zmin = WIN.[P].zmin[p]; }
  1946.   if (WIN.[P].zmax[p] != 1j) { zmax = WIN.[P].zmax[p]; }
  1947.   
  1948.   #
  1949.   # Check for potential errors
  1950.   #
  1951.   
  1952.   if (xmin == xmax) 
  1953.   { 
  1954.     # As good a guess as any
  1955.     xmin = xmin - 1;
  1956.     xmax = xmax + 1;
  1957.   }
  1958.   
  1959.   if (ymin == ymax) 
  1960.   { 
  1961.     # As good a guess as any
  1962.     ymin = ymin - 1;
  1963.     ymax = ymax + 1;
  1964.   }
  1965.   
  1966.   if (zmin == zmax) 
  1967.   { 
  1968.     # As good a guess as any
  1969.     zmin = zmin - 1;
  1970.     zmax = zmax + 1;
  1971.   }
  1972.   
  1973.   #
  1974.   # Finally, adjust if log-scales
  1975.   #
  1976.   
  1977.   if (find_char (WIN.[P].grid3x[p], "l"))
  1978.   {
  1979.     if (xmin <= 0 || xmax <= 0) { error ("plot: cannot plot log(x<=0)"); }
  1980.     xmin = log10 (xmin);
  1981.     xmax = log10 (xmax);
  1982.   }
  1983.   
  1984.   if (find_char (WIN.[P].grid3y[p], "l"))
  1985.   {
  1986.     if (ymin <= 0 || ymax <= 0) { error ("plot: cannot plot log(y<=0)"); }
  1987.     ymin = log10 (ymin);
  1988.     ymax = log10 (ymax);
  1989.   }
  1990.   
  1991.   if (find_char (WIN.[P].grid3z[p], "l"))
  1992.   {
  1993.     if (zmin <= 0 || zmax <= 0) { error ("plot: cannot plot log(z<=0)"); }
  1994.     zmin = log10 (zmin);
  1995.     zmax = log10 (zmax);
  1996.   }
  1997.   
  1998.   return 1;
  1999. };
  2000.  
  2001. ##############################################################################
  2002. #
  2003. # Find the X and Y scales for a list of matrices
  2004. #
  2005.  
  2006. list_scales = function ( data, key, p, Xmin, Xmax, Ymin, Ymax )
  2007. {
  2008.   once = 1;
  2009.   
  2010.   for (i in members (data))
  2011.   {
  2012.     M = real (data.[i]);
  2013.     if (class (M) != "num") { continue; }
  2014.     
  2015.     if (abs (key) > M.nc)
  2016.     {
  2017.       error_1 ("plot: KEY argument > M.nc");
  2018.     }
  2019.     
  2020.     #
  2021.     # 1st check for un-plottable data
  2022.     #
  2023.     
  2024.     if (any (any (isinf (M))))
  2025.     { error ("plot: cannot plot infs"); }
  2026.     if (any (any (isnan (M))))
  2027.     { error ("plot: cannot plot NaNs"); }
  2028.     
  2029.     k = find ((1:M.nc) != abs (key));
  2030.     if (key > 0)
  2031.     {
  2032.       if (M.nc != 1)
  2033.       {
  2034.     x_scales ( real(M)[;key], p, xmin, xmax );
  2035.     y_scales ( real(M)[;k],   p, ymin, ymax );
  2036.       else
  2037.     x_scales ( (1:M.nr)', p, xmin, xmax );
  2038.     y_scales ( real(M),   p, ymin, ymax );
  2039.       }
  2040.     else if (key < 0) {
  2041.       x_scales ( real(M)[;k],        p, xmin, xmax );
  2042.       y_scales ( real(M)[;abs(key)], p, ymin, ymax );
  2043.     else
  2044.       x_scales ( (1:M.nr)', p, xmin, xmax );
  2045.       y_scales ( real(M),   p, ymin, ymax );
  2046.     } }
  2047.  
  2048.     if (once) 
  2049.     { 
  2050.       Xmin = xmin; Xmax = xmax; Ymin = ymin; Ymax = ymax; 
  2051.       once = 0; 
  2052.     }
  2053.     if (xmin < Xmin) { Xmin = xmin; }
  2054.     if (xmax > Xmax) { Xmax = xmax; }
  2055.     if (ymin < Ymin) { Ymin = ymin; }
  2056.     if (ymax > Ymax) { Ymax = ymax; }
  2057.   }
  2058.   
  2059.   return 1;
  2060. };
  2061.  
  2062. ##############################################################################
  2063. #
  2064. # Find the maximum number of elements in a bin for a single 
  2065. # column matrix.
  2066. #
  2067.  
  2068. hist_scales = function ( data, nbin )
  2069. {
  2070.   dmin = min (real (data));
  2071.   dmax = max (real (data));
  2072.   dbin = linspace (dmin, dmax, nbin+1);
  2073.   binval = zeros (nbin, 1);
  2074.   
  2075.   for (i in 1:nbin)
  2076.   {
  2077.     binval[i] = length (find (data >= dbin[i] && data < dbin[i+1]));
  2078.   }
  2079.   
  2080.   return max (binval);
  2081. };
  2082.  
  2083. ##############################################################################
  2084. #
  2085. # Plot the columns of a matrix (core function)
  2086. #
  2087. # Notes: This is the core function for plotting a matrix. If the
  2088. # matrix is a single column, then the matrix elements are plotted
  2089. # versus the row numbers. If it is a multi-column matrix, then
  2090. # columns 2:N are plotted versus column 1.
  2091. #
  2092. # p, K, k and l are indices for plot features.
  2093. # p: the current plot index (the plot #)
  2094. # K: usually 0. This index is used to start of the line style and
  2095. # color index (k = color index, l = line-style index). This is mostly
  2096. # used by plot_list, which may call plot_matrix repeatedly.
  2097. # k: the line color index. This value determines the line color used
  2098. # for each column of data. If K = 0, then k goes like 2:14, then
  2099. # flops back to 1:14.
  2100. # l: the line style inex. This value determines the line style used
  2101. # for each column of data - not the line-type (points, or lines). If
  2102. # K = 0, then l goes like 2:8, then flops back to 1:8.
  2103. #
  2104. ##############################################################################
  2105.  
  2106. plot_matrix = function ( M, key, p, K, xmin, xmax, ymin, ymax, v )
  2107. {
  2108.   np = M.nr;
  2109.   
  2110.   if (M.nc == 1)
  2111.   {
  2112.     x = 1:M.nr;
  2113.     y = real (M);
  2114.     k = mod (1+K, 14) + 1;
  2115.     l = mod (1+K, 8);
  2116.     
  2117.     if (find_char (WIN.[P].gridx[p], "l"))
  2118.     { x = log10 (x); }
  2119.     if (find_char (WIN.[P].gridy[p], "l"))
  2120.     { y = log10 (y); }
  2121.     
  2122.     _plcol (WIN.[P].color[p;k]);
  2123.     _pllsty (WIN.[P].lstyle[p;l]);
  2124.     
  2125.     if (get_style (WIN.[P].style.[p], k-1) == "line") 
  2126.     {
  2127.       _plline (M.nr, x, y);
  2128.     else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2129.       _plpoin (M.nr, x, y, WIN.[P].pstyle[p]+k);
  2130.     else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2131.       _plline (M.nr, x, y);
  2132.       _plpoin (M.nr, x, y, WIN.[P].pstyle[p]+k);
  2133.     else {
  2134.       _plline (M.nr, x, y);
  2135.     }}}}
  2136.  
  2137.     #
  2138.     # Now do the legend 
  2139.     #
  2140.     
  2141.     if (!any (any (WIN.[P].desc.[p] == 1j)))
  2142.     {
  2143.       # Use the default if necessary
  2144.       if (WIN.[P].desc.[p][1] == "default") 
  2145.       {
  2146.     desc = "c1";
  2147.       else if (WIN.[P].desc.[p].n >= k-1) {
  2148.     desc = WIN.[P].desc.[p][k-1];
  2149.       else
  2150.     # Not sure what to do, user has messed up.
  2151.     desc = "";
  2152.       } }
  2153.              
  2154.       v = v - (ymax-ymin)/11;
  2155.       xl = (xmax-xmin)*[10.5/12, 11/12, 11.5/12]' + xmin;
  2156.       yl = [v, v, v]' + ymin;
  2157.       
  2158.       if (get_style (WIN.[P].style.[p], k-1) == "line") 
  2159.       {
  2160.     _plline (3, xl, yl);
  2161.       else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2162.     _plpoin (3, xl, yl, WIN.[P].pstyle[p]+k);
  2163.       else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2164.     _plline (3, xl, yl);
  2165.     _plpoin (3, xl, yl, WIN.[P].pstyle[p]+k);
  2166.       } } }
  2167.  
  2168.       plptex(desc, xl[1]-(xmax-xmin)/25, yl[3], , , 1);
  2169.       
  2170.     }
  2171.  
  2172.   else
  2173.  
  2174.     #
  2175.     # Check for large column dimension
  2176.     #
  2177.     
  2178.     if (M.nc > 3*M.nr)
  2179.     {
  2180.       printf (" Plot %i columns and %i rows, are you sure (y/n) ? "...
  2181.                , M.nc, M.nr);
  2182.       ans = getline ("stdin");
  2183.       if (ans.[1] != "y") { return -1; }
  2184.     }
  2185.     
  2186.     ki = find ((1:M.nc) != abs (key));
  2187.     for (i in ki)
  2188.     {
  2189.       if (key > 0)
  2190.       {
  2191.     x = real (M[;key]);
  2192.     y = real (M[;i]);
  2193.       else if (key < 0) {
  2194.     x = real (M[;i]);
  2195.     y = real (M[;abs(key)]);
  2196.       else
  2197.     x = (1:M.nr)';
  2198.     y = real (M[;i]);
  2199.       } }
  2200.  
  2201.       # Check for log scales, adjust if necessary
  2202.       if (find_char (WIN.[P].gridx[p], "l"))
  2203.       { x = log10 (x); }
  2204.       if (find_char (WIN.[P].gridy[p], "l"))
  2205.       { y = log10 (y); }
  2206.       
  2207.       k = mod (i-1 + K, 14) + 1;
  2208.       l = mod (8 + i-2 + K, 8) + 1;
  2209.       
  2210.       _plcol (WIN.[P].color[p;k]);
  2211.       _pllsty (WIN.[P].lstyle[p;l]);
  2212.       
  2213.       if (get_style (WIN.[P].style.[p], k-1) == "line") 
  2214.       {
  2215.     _plline (np, x, y);
  2216.       else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2217.     _plpoin (np, x, y, WIN.[P].pstyle[p]+k);
  2218.       else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2219.     _plline (np, x, y);
  2220.     _plpoin (np, x, y, WIN.[P].pstyle[p]+k);
  2221.       else {
  2222.     _plline (np, x, y);
  2223.       }}}}
  2224.  
  2225.       #
  2226.       # Now do the legend 
  2227.       #
  2228.       
  2229.       if (!any (any (WIN.[P].desc.[p] == 1j)))
  2230.       {
  2231.     # Use the default if necessary
  2232.     if (WIN.[P].desc.[p][1] == "default") 
  2233.     {
  2234.       desc = "c" + num2str (i);
  2235.         else if (WIN.[P].desc.[p].n >= k-1) {
  2236.       desc = WIN.[P].desc.[p][k-1];
  2237.         else
  2238.       # Not sure what to do, user has messed up.
  2239.       desc = "";
  2240.         }}
  2241.              
  2242.     v = v - (ymax-ymin)/11;
  2243.     xl = (xmax-xmin)*[10.5/12, 11/12, 11.5/12]' + xmin;
  2244.     yl = [v, v, v]' + ymin;
  2245.     
  2246.     if (get_style (WIN.[P].style.[p], k-1) == "line") 
  2247.     {
  2248.       _plline (3, xl, yl);
  2249.         else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2250.       _plpoin (3, xl, yl, WIN.[P].pstyle[p]+k);
  2251.         else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2252.       _plline (3, xl, yl);
  2253.       _plpoin (3, xl, yl, WIN.[P].pstyle[p]+k);
  2254.         }}}
  2255.     
  2256.     plptex(desc, xl[1]-(xmax-xmin)/25, yl[3], , , 1);
  2257.     
  2258.       }
  2259.     }
  2260.   }
  2261.   
  2262.   return k-1;
  2263. };
  2264.  
  2265. ##############################################################################
  2266. #
  2267. # Plot all of the matrices in a list on the same plot
  2268. #
  2269.  
  2270. plot_list = function ( L, key, p, xmin, xmax, ymin, ymax )
  2271. {
  2272.   k = 0;
  2273.   v = ymax - ymin;
  2274.   
  2275.   #
  2276.   # Sort out the list members
  2277.   #
  2278.  
  2279.   sl = list_sort (L);
  2280.  
  2281.   # Plot the list members with numeric labels 1st.
  2282.   if (exist (sl.num))
  2283.   {
  2284.     for (i in sl.num)
  2285.     {
  2286.       M = L.[i];
  2287.       if (class (M) != "num") { continue; }
  2288.       if ((k = plot_matrix (M, key, p, k, xmin, xmax, ymin, ymax, v)) < 0) 
  2289.       { 
  2290.     return k; 
  2291.       }
  2292.     }
  2293.   }
  2294.  
  2295.   # Now plot the list members with string labels.
  2296.   if (exist (sl.char))
  2297.   {
  2298.     for (i in sl.char)
  2299.     {
  2300.       M = L.[i];
  2301.       if (class (M) != "num") { continue; }
  2302.       if ((k = plot_matrix (M, key, p, k, xmin, xmax, ymin, ymax, v)) < 0) 
  2303.       { 
  2304.     return k; 
  2305.       }
  2306.     }
  2307.   }
  2308.   return 1;
  2309. };
  2310.  
  2311. ##############################################################################
  2312. #
  2313. # Check the elements of LIST.
  2314. # LIST must contain elements `x', `y',
  2315. # and `z'
  2316. #
  2317.  
  2318. check_3d_list = function ( LIST )
  2319. {
  2320.   #
  2321.   # Check existence and types
  2322.   #
  2323.   
  2324.   if (class (LIST) != "list") {
  2325.     error ("plot3: argument must be a list");
  2326.   }
  2327.   if (!exist (LIST.x)) {
  2328.     error ("plot3: arg must contain `x' member");
  2329.   else if (class (LIST.x) != "num") {
  2330.     error ("plot3: x must be numeric");
  2331.   } }
  2332.   if (!exist (LIST.y)) {
  2333.     error ("plot3: arg must contain `y' member");
  2334.   else if (class (LIST.y) != "num") {
  2335.     error ("plot3: y must be numeric"); 
  2336.   } }
  2337.   if (!exist (LIST.z)) {
  2338.     error ("plot3: arg must contain `z' member");
  2339.   else if (class (LIST.z) != "num") {
  2340.     error ("plot3: z must be numeric");
  2341.   } }
  2342.  
  2343.   #
  2344.   # Check sizes
  2345.   #
  2346.   
  2347.   if (LIST.x.n != LIST.z.nr) 
  2348.   {
  2349.     error ("plot3: x.n != z.nr");
  2350.   }
  2351.   
  2352.   if (LIST.y.n != LIST.z.nc) 
  2353.   {
  2354.     error ("plot3: y.n != z.nc");
  2355.   }
  2356.   
  2357. };
  2358.  
  2359. ##############################################################################
  2360. #
  2361. # A special type of histogram plot.
  2362. #
  2363.  
  2364. plhistx = function ( M , nbin )
  2365. {
  2366.   check_plot_object ();
  2367.   
  2368.   if (!exist (nbin)) { nbin = 10; }
  2369.   
  2370.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  2371.   np = M.nr;
  2372.   
  2373.   # Compute max/min values of data
  2374.   
  2375.   ymin = min (min(real (M)));
  2376.   ymax = max (max(real (M)));
  2377.   
  2378.   #
  2379.   # Check computed scale limits against user's
  2380.   #
  2381.   
  2382.   if (WIN.[P].ymin[p] != 1j) { ymin = WIN.[P].ymin[p]; }
  2383.   if (WIN.[P].ymax[p] != 1j) { ymax = WIN.[P].ymax[p]; }
  2384.   
  2385.   _plgra ();
  2386.   _plcol (15);
  2387.   _plfont (WIN.[P].font[p]);
  2388.   _plwid (WIN.[P].width[p]);
  2389.   
  2390.   dbin = (linspace (ymin, ymax, nbin+1))';
  2391.   for (j in 1:M.nc) 
  2392.   {
  2393.     // counting
  2394.     for (i in 1:nbin) 
  2395.     {
  2396.       binval[i;j] = length (find (M[;j] >= dbin[i] && M[;j] < dbin[i+1]));
  2397.     }
  2398.   }
  2399.   
  2400.   if (!subplot_f) {
  2401.     _pladv (0);        # Advance 1 subplot
  2402.   else
  2403.     subplot_f = 0;     # The user has set the subplot
  2404.   }
  2405.  
  2406.   if (WIN.[P].aspect[p] != 0)
  2407.   {
  2408.     _plvasp (WIN.[P].aspect[p]);
  2409.   else
  2410.     _plvsta ();
  2411.   }
  2412.  
  2413.   xmin = 0;
  2414.   xmax =  max(max(binval));
  2415.   _plwind (ymin, ymax, xmin, xmax);
  2416.   _plbox (WIN.[P].gridx[p], 0, 0, WIN.[P].gridy[p], 0, 0);
  2417.  
  2418.   #
  2419.   # Now reorganize dbin and binval so we plot points in
  2420.   # the middle of the bins.
  2421.   #
  2422.  
  2423.   delbin = abs(dbin[2] - dbin[1]);
  2424.   dbin = (linspace (ymin+delbin/2, ymax-delbin/2, nbin))';
  2425.   dbin = [ymin ; dbin ; ymax];
  2426.   binval = [zeros(1,binval.nc); binval; zeros(1,binval.nc)];
  2427.  
  2428.   v = xmax;
  2429.   for (i in 1:M.nc)
  2430.   {
  2431.     k = mod (i, 14) + 1;
  2432.     l = mod (i,  8) + 1;
  2433.     _plcol (WIN.[P].color[p;k]);
  2434.     _pllsty (WIN.[P].lstyle[p;l]);
  2435.     
  2436.     if      (get_style (WIN.[P].style.[p], k-1) == "line") 
  2437.     {
  2438.       _plline (nbin+2, dbin, binval[;i]);
  2439.     else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2440.       _plpoin (nbin+2, dbin, binval[;i], WIN.[P].pstyle[p]+k);
  2441.     else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2442.       _plline (nbin+2, dbin, binval[;i]);
  2443.       _plpoin (nbin+2, dbin, binval[;i], WIN.[P].pstyle[p]+k);     
  2444.     } } }
  2445.  
  2446.     // write legend around upper-right corner.
  2447.     // it is better to have user to choose location for legend.
  2448.     
  2449.     if (!any (any (WIN.[P].desc.[p] == 1j)))
  2450.     {
  2451.       # Use the default if necessary
  2452.       if (WIN.[P].desc.[p][1] == "default") 
  2453.       {
  2454.     desc = "c"+num2str(i);
  2455.       else if (WIN.[P].desc.[p].n >= i) {
  2456.       desc = WIN.[P].desc.[p][i];
  2457.       else
  2458.     # Not sure what to do, user has messed up.
  2459.     desc = "";
  2460.       } }
  2461.  
  2462.       v = v - (xmax)/11;
  2463.       xt = (ymax-ymin)*[10.5/12, 11/12, 11.5/12]' + ymin;
  2464.       yt = [v, v, v]';
  2465.       
  2466.       if      (get_style (WIN.[P].style.[p], k-1) == "line") 
  2467.       {
  2468.     _plline (3, xt, yt);
  2469.       else if (get_style (WIN.[P].style.[p], k-1) == "point") {
  2470.     _plpoin (3, xt, yt, WIN.[P].pstyle[p]+k);
  2471.       else if (get_style (WIN.[P].style.[p], k-1) == "line-point") {
  2472.     _plline (3, xt, yt);
  2473.     _plpoin (3, xt, yt, WIN.[P].pstyle[p]+k);
  2474.       } } }
  2475.  
  2476.       plptex(desc, xt[1]-(ymax-ymin)/25, yt[3], , , 1);
  2477.     }
  2478.   }
  2479.   
  2480.   _plcol (15);
  2481.   _pllab (WIN.[P].xlabel[p], WIN.[P].ylabel[p], WIN.[P].title[p]);  
  2482.   _plflush ();
  2483.   _pltext ();
  2484.   
  2485.   #
  2486.   # Increment the plot no. so that next time
  2487.   # we use the correct settings.
  2488.   #
  2489.   
  2490.   WIN.[P].subplot = WIN.[P].subplot + 1;
  2491.   
  2492.   return 1;
  2493. };
  2494.  
  2495.  
  2496. ##############################################################################
  2497. #
  2498. # Create a legend in the current plot window
  2499. #
  2500. # if pobj.desc.[p] = inf()        no legend
  2501. # if pobj.desc.[p] = "default"        default ("c1", "c2", ...)
  2502. # if pobj.desc.[p] = "string"        use "string" as description
  2503. #
  2504.  
  2505. #
  2506. # Set the current plot legend string
  2507. #
  2508.  
  2509. plegend = function ( LEGEND )
  2510. {
  2511.   check_plot_object ();
  2512.   
  2513.   p = mod (WIN.[P].subplot, WIN.[P].nplot) + 1;    # The current index
  2514.   
  2515.   if (!exist (LEGEND)) 
  2516.   {
  2517.     WIN.[P].desc.[p] = 1j;
  2518.     return P;
  2519.   }
  2520.   
  2521.   if (class (LEGEND) == "string")
  2522.   {
  2523.     WIN.[P].desc.[p] = LEGEND;
  2524.   }
  2525.   
  2526.   return P;
  2527. };
  2528.  
  2529. set3d = function (bx, by, h)
  2530. {
  2531.   if (!exist (bx)) { basex = 2; else basex = bx; }
  2532.   if (!exist (by)) { basey = 2; else basey = by; }
  2533.   if (!exist (h)) { height = 4; else height = h; }
  2534. };
  2535.